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We found an analytical solution for the vortex structure in a rapidly rotating trapped Bose- 
Einstein condensate in the lowest Landau level approximation. This solution is exact in the limit of 
a large number of vortices and is obtained for the case of anisotropic harmonic potential. For the 
case of symmetric harmonic trap when the rotation frequency is equal to the trapping frequency, 
the solution coincides with the Abrikosov triangle vortex lattice in type-II superconductors. In a 
general case the coarse grained density is found to be close to the Thomas-Fermi profile, except the 
vicinity of edges of a condensate cloud. 
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Bose- Einstein condensates (BEC) are a new state of matter where various aspects of macroscopic quantum physics 
can be studied. The experimental studying of BEC in ultra-cold rotating atomic gases shows a wide variety of new 
features in the physics of quantized vortices and vortex arrays [Il-Q which can be not accessible in other quantum 
systems containing vortices, such as superfluid helium or type-II superconductors. In a harmonically trapped conden- 
sate rotating at a frequency close to the trap frequency, vortices form a triangular Abrikosov lattice and the coarse 
grained density approaches an inverted parabola [^-l7| . At very fast rotation, when the number of vortices becomes 
close to the number of atoms, the states are strongly correlated and the vortex lattice is expected to melt. 

We propose below the solution for the vortex structure of rotating BEC in anisotropic harmonic trapping. The found 
function satisfiers the projected to the lowest Landau level (LLL) Gross-Pitaevskii (CP) equation with extremely high 
accuracy for a large number of vortices. 

We consider a system of bosonic atoms strongly confined in the z direction by an external trapping potential such 
that the bosons become essentially two-dimensional. The bosons are confined in the plane by a anisotropic harmonic 
trapping potential ^(r) = m(a;^a;^ -I- a;^y^)/2, with ujy < iOx, for a definiteness, and the trap is rotating around the 
z axis with frequency fi. We assume that all particles are in the same macroscopic quantum state described by the 
wave function '(/'(r)- In the rotating frame the Gross-Pitaevskii equation for 'ip{v) reads: 
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where p is the momentum operator, m is the particle mass, 5 > is an effective 2D coupling constant, is the 
operator of the orbital angular momentum, fi is the chemical potential, and ip is normalized to the total number of 
particles N. 

As it is well known, the single-particle Hamiltonian Hq for rotating neutral atoms is equivalent to the Hamiltonian 
of a charged particle in a uniform magnetic field B — 2mfl along the z axis. The half the cyclotron frequency is equal 
to Uc — ri, and the vector-potential in the symmetric gauge is A = B x r/2 = mO x r. 
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Below the critical rotation frequency, fl < uj^, LOy, the presence of the residual confining potential lifts the degeneracy 



of the Landau levels. Provided that i 



-Q fl and interactions are mach smaller than the cyclotron frequency 



{ng ^ 2fl, where n is the two-dimensional particle density), we can restrict our consideration to the lowest Landau 
level. 

The Gross-Pitaevskii equation projected onto the LLL of an asymmetric harmonic trap can be written in terms of 
dimensionless variables 5, y as 
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with ^2 ^ a;2 - ujj^ (1)2 ^ a;2 - ujj, ujf 
have a form 



^ u):^±uy g^j^j^ dxdy = (wj'^Wj /u)c)dxdy. The LLL wave functions 



^{x,y) ^ /(C)e-5"t ^s^-2(^^ 
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with an analytic function /(C), and a = \p'^{ijj^ — ujc), b = h{(^t + ^c)l 
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The anisotropy parameter tanhp is equal to zero for the cyhnder symmetric trap: = uiy > fl. The opposite 
quasi-one-dimensional case with a narrow channel geometry is achieved at tanhi/ — >■ 1, when the rotation frequency 
becomes equal to the lowest trapping frequency: fl ^ ujy < lj^- Both limited cases were considered in details in 

Note that for the case oi ujx = i^y = ^ the free energy J dx dy[ip* Hoip + g\ip\'^ /2 — has the same form as 

the Ginzburg-Landau functional of a superconductor in the magnetic field near a phase transition. The important 
difference and simplification in our case is the constant "magnetic field" 2m^ (infinite penetration depth). Therefore, 
instead of three equations obtained by variation over ip, B, plus the boundary condition, we have one equation with the 
normalization condition J dxdy\ip'^\ = N. Since the penetration depth of the "magnetic field" is infinity, properties 
of our system will be similar to properties of type-II superconductors. Therefore it is reasonable to expect that 
superconductors in the magnetic field and rotating BEC will have similar structure of vortex lattice. Indeed, we can 
see that the well known approximate solution for Abrikosov vortex lattice Q was built by use of the lowest Landau 
level (LLL) wave functions, as well as the solution for the BEC. Moreover, an approximate solution for Abrikosov 
lattice in type-II superconductor becomes the exact solution for the considered model of rotating BEC in the LLL 
approximation. 

For the general case of asymmetric harmonic potential the ground state contains an ordered vortex structure in the 
parameters region Ng/P /{h{ujx^y — fl)) ^ 1, where / = {h/mfiy^^ is the effective magnetic length. The number of 
vortices will increase with the increase of this ratio. When the number of particles will become of the order of the 
number of vortices, the vortex lattice will melt through the phase transition to a strongly correlated state. In this 
region the Gross-Pitaevskii mean field equation is not applicable. 

We will find the solution for the vortex structure as a special deformation of the exact solution for spatially 
homogeneous vortex lattice for the case of a cylinder symmetric potential (lo^ = ujy = lj) at the critical value of the 
frequency rotation Q. — io: 

/o(C) = ^^^'^i {</bur) exp(^CV2«c), (7) 
where S is the surface area, and di is the Jacobi theta-function given by 

oo 

di{r^,T)^- y (-l)"exp{^7rT(r^ + 1/2)2 + 2i7;(n + 1/2)}. 

n— — oo 

Real (m) and imaginary {iv) parts of the quasiperiod t = u + iv will be fixed below. The Jacobi theta-functions are 
analytic in the complex plane and have zeros at the points "q = mr + rmrT, where n, m are integers. The function 
/o(C) has zeros at the lattice sites nbi + mb2, with &2 — biT. These points correspond to the vortex locations. The 
parameter hi can be chosen real so that the area of the unit cell is Vc = h\v. Using the property of the t?-function: 



K/&i,r)| =G(i,y)exp 



with a periodic oscillating G: < G{x,y) < §2 (0,t), one can obtain that the envelope |/o(C)l has a polar symmetric 
form 



l/o(C)l^exp 
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The function fo{z) with fixed elementary cell area = tt is an exact solution of Eq. ^ for the case of a critical rotation 
cja; = cjy = fl. The normalization coefficient in Eq. ([7]) is chosen such that the function 5* ~ [fo{z)/l] exp(— zz/2) is 
normalized to unity. The function ^ has a constant envelope and describes a periodic vortex structure. The minimum 
energy is obtained for the triangular lattice, where r = exp27ri/3, v = and 6^ — 2n/^/3. The chemical potential 

is then given by = with a = (31/^/2) ^^^^(-1)"^ exp{-7r2(fe2 + c^yAbj} = 1.1596, and b ^ 2m, c = 2p being 
even integers. 

In the general case ( fl < u!x,^y), we will find a solution of equation ([3]) in the form 

/(C) = (2v)i/4 J2 (-l)"5(a)g'^'e*v^"«+V(i-*'^"h^\ (10) 

n— — oo 

where a = 2n -I- 1 are odd integers, q = exp[i7rr/4], and g{a) is a differential operator acting on a. Substituting the 
trial function (ITU|) into equation for a triangular-like lattice (the lattice that becomes exactly triangular for lu = n) 
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we obtain: 



X cxp 
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Here we introduced the operators A, A' which a creation (annihilation) operators for a corresponding harmonic 
osciUator with usual commutation relations [A. A^] = 1: 
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firfi da 



where v — vj, with 7 = (1 

For large /3 = Ng/{P{uj^ 
for g{a), which describes the vortex structure with a high accuracy, has the form 



tanh v)/{l — tanh v) . 

)) 3> 1 and /i* = (/i — uj^)/{u}^ — uj^) — 1/(1 + tanh i^) 3> 1 an approximate solution 
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where Q is the Heaviside step function, R ~ y/Jj^, and a = 1.1596. We will see below that for the symmetric case 
R is the radius of the condensate cloud in units of /. Equation is obtained taking into account that the leading 
contribution to the sum over b and c in Eq. (1111) comes from small values of m and p, since already the contributions 
of terms with |m| > 2 or \p\ > 2 are exponentially small. Provided that the dependence g{a) is smooth, which is the 
case for large R, we may consider large a and omit b and c in the arguments of the g operators in Eq. (jlll) . This 
immediately gives Eq. (1121) . Mathematically, the high accuracy of the solution is based on a known fast convergency 
of the series for ??-function due to the multiplier gt-^n+i) _ 



Substituting the solution (fT^ into Eq. 
oscillator we obtain after some algebra: 



TU)) and expanding (|10p into series in known eigenfunctions of harmonic 



j-^^^ ^ (2«)^/V(l + tanh;/) g ^^_-^^[„(„-i)/2] 



VR^-k 



(iVtanhz^)*^ 



ri=— 00 k=0 



2*=fc! 



f c 



\ vsinh2l' 



(13) 




where Hk{w) are Hermite polynomials, a = 2n + 1. The solution (|13p is simplified in the symmetric case, where 
— )■ 0, -fffc(C/Vsinh2i^) — >■ 2^/"^^^ /v^^"^ ^ and in the one-dimensional case where v ^ 00 and operators A, ^4^' becomes 
numbers: A = A^ tx a. 

From the condition that the function [/(C)/'] exp(— CC) is normalized to unity we find a relation R = (2afi^ jii)^^'^ ■ 
For the symmetric potential (7 = 1) this result is in agreement with Refs. [3, S [9l-[ll| . As mentioned above the 
solution JlSp includes limiting cases of cylinder and narrow channel geometries, considered in details in 7]. Numerical 
results [7| for these cases demonstrated excellent coincidence with the analytical solution. The structure of the vortex 
lattice for i? = 7, tanhi/ = 1/4 is shown in Fig. [TJ Averaging the density over the oscillations, that is averaging the 
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FIG. 1. (Color online) Condensate wave-function \ip{x,y)f for R — 7, tanhi/ — 1/4. Coordinates x and y are given in units of 
/. 

density j^'p over a distance scale much larger than I, gives the coarse grained density: 
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r(i?2 + 1) ^ 



where ^(^1,^2) = 2 tanh(z/)(— i(/-\/sinh(2i/) + M(i(^/-y/sinh(2i^) +t2)- In the limit i/ ^ when X — >• we obtain 
exactly the expression for the symmetric case [7|. The solution for a narrow channel case is obtained in the limit 
I' — CXI taking into account the divergence of i?: i? (x (1 — tanhi^)"^/^. 

The expression can be simplified in the case of large R. Inside the condensate cloud we obtain the density 
profile which is close to the Tomas-Fermi inverted paraboloid: 



n = ncg(O) 



tanh v 



i?2(l-tanh^i/ 
+ i^(l - tanhi^)^ + 2/^(1 + tanh z/)^] 



(15) 



with ncg(O) ~ 2N / {-kR ). The condensate is localized at the region inside the ellipse x jar + y /b = 1 with semiaxes 
a>b: 



R_ 



Outside the ellipse the density is negligible small and has the asymptotic form 

TT-cg ocexp[i?^ ln(|CP(l -tanhi/)^) 
-i^(l - tanhi^) - y^(l +tanhi^)]. 



(16) 
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Outside the condensate cloud the density drops exponentially at distance of the order of magnetic length or intervortex 
spacing. In the direction y, for example, in the interval (1 <C y — & ^ b^^^) we obtain 



"cg(0) 



-2(y-bf 



(17) 
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FIG. 2. (Color online) The coarse grained density Ucgix, 0) in units of n^D versus for R — 5, taiii/ = 1/2. The dashed curve 
shows the Thomas-Fermi inverted-parabola shape, and x is given in units of /. 
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FIG. 3. (Color online) The coarse grained density ncg(0, y) in units of 7120 versus y, for R ■ 
shows the Thomas-Fermi inverted-parabola shape, and y is given in units of I. 



5, tanz/ = 1/2. The dashed curve 



At the Thomas-Fermi border we have ricg ~ ncg(0)/-\/27r|C|2. Our results show that inside the condensate cloud 
the density profile has the Thomas- Fermi inverted-paraboloid shape, except the region close to the edge: \Sz\ < I, 
where \Sz\ is the distance to the border. We see in Figs. [U [3] that for R = 5, tanz/ = 1/2 the Thomas-Fermi 
formula works well already for \SC\ ^ 1- The parabolic coarse grained density profile for the asymmetric potential was 
considered analytically and numerically in |12h14| . Deviations from the Thomas- Fermi density profile of ncg{r) have 
been studied in Ref. [5| for the symmetric potential by using the variational numeric methods. Here we presented an 
analytic solution for the vortex lattice in a rapidly rotating BEC in an asymmetric harmonic trap. The found coarse 
grained density profile is close to inverted paraboloid form. The solution is asymptotically exact in the limit of a large 
number of vertices. 
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